//
// Created on 2025/5/29/029.
//
// Node APIs are not fully supported. To solve the compilation error of the interface cannot be found,
// please include "napi/native_api.h".


//#include "YIN.h"
#include "napi/native_api.h"
#include <cstdint>
#include <math.h>
#include <hilog/log.h>

const double PI = 3.14159265358979323846;

typedef struct {
    float *buffer;          // 原始音频数据
    float *diff;            // 差分函数结果
    float *cmndf;           // 累加归一化差分函数
    size_t sample_size;   // 帧大小（建议2048）
    float threshold;        // 阈值（如 0.15）
    float sample_rate;      // 采样率（如 16000）

    float min_freq;       // 可检测的最小频率
    float max_freq;       // 可检测的最大频率
    uint32_t min_tau;     // 可检测的最小周期
    uint32_t max_tau;     // 可检测的最大周期
} YIN_Instance;

static YIN_Instance *yin_instance;

//static void arm_power_f32(float *a,uint32_t len,float *result)
//{
//    
//}

//static void arm_dot_prod_f32(float *a,float *b,uint32_t len,float *result)
//{
//    
//}

typedef float float32_t;
#define ARM_MATH_LOOPUNROLL

void arm_power_f32(
  const float32_t * pSrc,
        uint32_t blockSize,
        float32_t * pResult)
{
        uint32_t blkCnt;                               /* Loop counter */
        float32_t sum = 0.0f;                          /* Temporary result storage */
        float32_t in;                                  /* Temporary variable to store input value */

  /* Loop unrolling: Compute 4 outputs at a time */
  blkCnt = blockSize >> 2U;

  while (blkCnt > 0U)
  {
    /* C = A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1] */

    /* Compute Power and store result in a temporary variable, sum. */
    in = *pSrc++;
    sum += in * in;

    in = *pSrc++;
    sum += in * in;

    in = *pSrc++;
    sum += in * in;

    in = *pSrc++;
    sum += in * in;

    /* Decrement loop counter */
    blkCnt--;
  }

  /* Loop unrolling: Compute remaining outputs */
  blkCnt = blockSize % 0x4U;


  while (blkCnt > 0U)
  {
    /* C = A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1] */

    /* Compute Power and store result in a temporary variable, sum. */
    in = *pSrc++;
    sum += in * in;

    /* Decrement loop counter */
    blkCnt--;
  }

  /* Store result to destination */
  *pResult = sum;
}

void arm_dot_prod_f32(
  const float32_t * pSrcA,
  const float32_t * pSrcB,
        uint32_t blockSize,
        float32_t * result)
{
        uint32_t blkCnt;                               /* Loop counter */
        float32_t sum = 0.0f;                          /* Temporary return variable */

  /* Loop unrolling: Compute 4 outputs at a time */
  blkCnt = blockSize >> 2U;

  /* First part of the processing with loop unrolling. Compute 4 outputs at a time.
   ** a second loop below computes the remaining 1 to 3 samples. */
  while (blkCnt > 0U)
  {
    /* C = A[0]* B[0] + A[1]* B[1] + A[2]* B[2] + .....+ A[blockSize-1]* B[blockSize-1] */

    /* Calculate dot product and store result in a temporary buffer. */
    sum += (*pSrcA++) * (*pSrcB++);

    sum += (*pSrcA++) * (*pSrcB++);

    sum += (*pSrcA++) * (*pSrcB++);

    sum += (*pSrcA++) * (*pSrcB++);

    /* Decrement loop counter */
    blkCnt--;
  }

  /* Loop unrolling: Compute remaining outputs */
  blkCnt = blockSize % 0x4U;

  while (blkCnt > 0U)
  {
    /* C = A[0]* B[0] + A[1]* B[1] + A[2]* B[2] + .....+ A[blockSize-1]* B[blockSize-1] */

    /* Calculate dot product and store result in a temporary buffer. */
    sum += (*pSrcA++) * (*pSrcB++);

    /* Decrement loop counter */
    blkCnt--;
  }

  /* Store result in destination buffer */
  *result = sum;
}



// 差分函数计算
void yin_difference(YIN_Instance *yin)
{
    uint32_t N = yin->sample_size;
    float *x = yin->buffer;
    float *d = yin->diff;

    d[0] = 0.0f;

    for (uint32_t tau = 1; tau < N / 2; tau++)
    {
        uint32_t length = N - tau;

        float energy1 = 0.0f;
        float energy2 = 0.0f;
        float cross = 0.0f;

        arm_power_f32(&x[0], length, &energy1);
        arm_power_f32(&x[tau], length, &energy2);
        arm_dot_prod_f32(&x[0], &x[tau], length, &cross);

        d[tau] = energy1 + energy2 - 2.0f * cross;
    }
}

void yin_difference_optimized(YIN_Instance *yin)
{
    uint32_t N = yin->sample_size;
    float *x = yin->buffer;
    float *d = yin->diff;
    float *cumulative_energy = yin->cmndf; // 复用缓冲区
    float energy1, energy2, cross;

    // 1. 预计算累积能量
    cumulative_energy[0] = x[0] * x[0];
    for (uint32_t i = 1; i < N; i++) {
        cumulative_energy[i] = cumulative_energy[i - 1] + x[i] * x[i];
    }

    float total_energy = cumulative_energy[N - 1];
    d[0] = 0.0f;

    // 2. 使用累积能量计算差分
    for (uint32_t tau = 1; tau < N / 2; tau++) {
        uint32_t length = N - tau;

        energy1 = cumulative_energy[length - 1];
        energy2 = total_energy - cumulative_energy[tau - 1];// x[tau] 到 x[N-1] 的能量
        cross = 0.0f;

        arm_dot_prod_f32(&x[0], &x[tau], length, &cross);

        d[tau] = energy1 + energy2 - 2.0f * cross;
    }
}


// 累加归一化
void yin_cumulative_mean_normalized_difference(YIN_Instance *yin)
{
    float running_sum = 0.0f;
    float *d = yin->diff;
    float *cmndf = yin->cmndf;
    uint32_t tau;

    for (tau = 2; tau < yin->sample_size / 2; tau++)
    {
        running_sum += d[tau];
        cmndf[tau] = d[tau] * tau / running_sum;
    }
    cmndf[0] = 1.0f;
    cmndf[1] = 1.0f;
}

// 寻找第一个小于阈值的极小值（使用 min_tau/max_tau 范围）
int yin_absolute_threshold(YIN_Instance *yin)
{
    uint32_t min_tau = yin->min_tau;
    uint32_t max_tau = yin->max_tau;

    if (max_tau >= yin->sample_size / 2)
    {
        max_tau = yin->sample_size / 2 - 1;
    }

    for (uint32_t tau = min_tau; tau <= max_tau; tau++)
    {
        if (yin->cmndf[tau] < yin->threshold)
        {
            // 向后找极小值
            while (tau + 1 <= max_tau && yin->cmndf[tau + 1] < yin->cmndf[tau])
            {
                tau++;
            }
            return tau;
        }
    }

    return -1; // 未找到周期
}

// 插值函数
float yin_parabolic_interpolation(float *cmndf, int tau)
{
    if (tau <= 1 || tau + 1 >= 1024)
        return (float)tau; // 边界检查

    float s0 = cmndf[tau - 1];
    float s1 = cmndf[tau];
    float s2 = cmndf[tau + 1];

    float denom = 2 * (s0 - 2 * s1 + s2);
    if (fabsf(denom) < 1e-6f)
        return (float)tau;

    float delta = (s0 - s2) / denom;
    return tau + delta;
}

// 应用于 yin->buffer 的 Hamming 窗
void apply_hamming_window(YIN_Instance *yin)
{
    uint32_t N = yin->sample_size;
    float *x = yin->buffer;

    for (uint32_t n = 0; n < N; n++)
    {
        float w = 0.54f - 0.46f * cos((2.0f * PI * n) / (N - 1));
        x[n] *= w;
    }
}

float yin_get_frequency(YIN_Instance *yin)
{
    if(yin == nullptr)
    {
        return -1.0f;
    }
    if(yin->buffer == nullptr)
    {
        return -1.0f;
    }
    // apply_hamming_window(yin);
    // yin_difference(yin);
    yin_difference_optimized(yin);
    yin_cumulative_mean_normalized_difference(yin);
    int tau = yin_absolute_threshold(yin);
    if (tau != -1)
    {
        // return yin->sample_rate / tau;
        float refined_tau = yin_parabolic_interpolation(yin->cmndf, tau);
        return yin->sample_rate / refined_tau;
    }
    return -1.0f; // 无法检测频率
}

void yin_init(float threshold, size_t frame_sample_size, float sample_rate, float min_freq, float max_freq)
{
    YIN_Instance *yin = new YIN_Instance;
    yin->sample_size = frame_sample_size;
//    yin->buffer = new float[yin->sample_size];
    yin->buffer = nullptr;
    yin->diff = new float[yin->sample_size];
    yin->cmndf = new float[yin->sample_size];
    yin->threshold = threshold;
    yin->sample_rate = sample_rate;
    yin->max_freq = max_freq;
    yin->min_freq = min_freq;
    yin->min_tau = (uint32_t)(yin->sample_rate / yin->max_freq);
    yin->max_tau = (uint32_t)(yin->sample_rate / yin->min_freq);
    yin_instance = yin;
}

void yin_free(YIN_Instance *yin)
{
    if(yin != nullptr)
    {
//        delete [] yin->buffer;
        delete [] yin->diff;
        delete [] yin->cmndf;
        delete yin;
    }
}

napi_value tuner_free(napi_env env, napi_callback_info info)
{
    yin_free(yin_instance);
    yin_instance = nullptr;
    
    return nullptr;
}

napi_value tuner_init(napi_env env, napi_callback_info info)
{
    size_t argc = 4;                // Number of parameters.
    napi_value args[4] = {nullptr}; // A napi array that stores parameters.

    if (napi_ok != napi_get_cb_info(env, info, &argc, args, nullptr, nullptr)) {
        OH_LOG_Print(LOG_APP, LOG_ERROR, 0xFF00, "ArraybufferType", "api_get_cb_info failed");
        return nullptr;
    }

    int32_t frameSampleSize;
    int32_t sampleRate;
    int32_t minFreq;
    int32_t maxFreq;
    if ((napi_ok != napi_get_value_int32(env, args[0], &frameSampleSize))
        || (napi_ok != napi_get_value_int32(env, args[1], &sampleRate))
        || (napi_ok != napi_get_value_int32(env, args[2], &minFreq)) 
        || (napi_ok != napi_get_value_int32(env, args[3], &maxFreq))) {
        OH_LOG_Print(LOG_APP, LOG_ERROR, 0xFF00, "ArraybufferType", "napi_get_arraybuffer_info failed");
        return nullptr;
    }
    
    yin_init(0.3,frameSampleSize,sampleRate,minFreq,maxFreq);

    return nullptr;
}

napi_value tuner_get_frequency(napi_env env, napi_callback_info info)
{
    size_t argc = 1;                // Number of parameters.
    napi_value args[1] = {nullptr}; // A napi array that stores parameters.

    if (napi_ok != napi_get_cb_info(env, info, &argc, args, nullptr, nullptr)) {
        OH_LOG_Print(LOG_APP, LOG_ERROR, 0xFF00, "ArraybufferType", "api_get_cb_info failed");
        return nullptr;
    }

    void *data = nullptr; // The pointer to the underlying data buffer used to get the arraybuffer.
    size_t length;        // The length of the underlying data buffer used to obtain the arraybuffer.

    napi_status res;
    res = napi_get_arraybuffer_info(env, args[0], &data, &length);
    if (napi_ok != res) {
        OH_LOG_Print(LOG_APP, LOG_ERROR, 0xFF00, "ArraybufferType", "napi_get_arraybuffer_info failed");
        return nullptr;
    }

    double freq = -1;
    size_t frameSize = length / sizeof(float);
    if(frameSize == yin_instance->sample_size)
    {
        yin_instance->buffer = reinterpret_cast<float *>(data);
        freq = yin_get_frequency(yin_instance);
    }

    napi_value frequency;
    if (napi_ok != napi_create_double(env, freq, &frequency)) {
        OH_LOG_Print(LOG_APP, LOG_ERROR, 0xFF00, "NumberType", "napi_create_int64 failed");
        return nullptr;
    }
    return frequency;
}